pacman::p_load(sf, tidyverse, tmap, dplyr,
spatstat, spdep,
lubridate, leaflet,
plotly, DT, viridis,
ggplot2, sfdep)Visual Analytics exercise on Armed conflicts - Initial Draft
4.1 The Task
In this take-home exercise, we are required to select one of the modules of our proposed Shiny application (Group Project) and complete the following tasks:
To evaluate and determine the necessary R packages needed for our Shiny application are supported in R CRAN,
To prepare and test that the specific R codes can run and returns the correct output as expected,
To determine the parameters and outputs that will be exposed on the Shiny applications,
To select the appropriate Shiny UI components for exposing the parameters determined above, and
To include a section called UI design for the different components of the UIs for the proposed design.
4.2 Getting Started
Our project will be using open-source data from the Armed Conflict Location & Event Data Project (ACLED).
Specifically, our project will be focusing on the visualisation of Armed conflicts in Myanmar, and I will be preparing the modules on Cluster & Outlier Analysis (LISA) and Hot/Cold zone analysis.
4.2.1 Loading R packages
The below R packages will be used in this exercise and for the Shiny application
4.2.2 Importing and loading the ACLED data
Country specific data from the Armed Conflict Location & Event Data Project (ACLED) can be downloaded at https://acleddata.com/data-export-tool/
Loading the ACLED data set for Myanmar.
ACLED_MMR <- read_csv("data/MMR.csv")4.2.3 Downloading and loading the shape files for country
Shape files were downloaded from the Myanmmar Information Management Unit (MIMU) website.
I chose this source over GADM and GeoBoundaries due to its updated administrative region information and map levels.
Note - Data Quality Issues
ACLED captures event data from national, sub-national and other media sources, and populates event locations based on the last known information.
However, some names of administrative areas were found to have changed; either disaggregated into new administrative areas or previously active but now defunct. Some administrative areas were also aggregated into higher administrative areas.
As part of our data cleaning and preparation process, I had to identify discrepancies in both admin1 & 2 (administrative levels) and re-name some administrative areas to sync with the downloaded shape files from MIMU.
4.3 Data Preparation and Cleaning
I will first load in the shape files at the admin 1 (region/state) and admin 2 (districts) levels. Most of our plots will be utilizing admin 2 levels.
4.3.1 Loading Admin 1 & 2 (administrative region/area) shape files
mmr_shp_mimu_1 <- st_read(dsn = "data/geospatial3",
layer = "mmr_polbnda2_adm1_250k_mimu_1")Reading layer `mmr_polbnda2_adm1_250k_mimu_1' from data source
`C:\imranmi\imran's data sc\R-ex\R-Ex5\data\geospatial3' using driver `ESRI Shapefile'
Simple feature collection with 18 features and 6 fields
Geometry type: MULTIPOLYGON
Dimension: XY
Bounding box: xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS: WGS 84mmr_shp_mimu_2 <- st_read(dsn = "data/geospatial3",
layer = "mmr_polbnda_adm2_250k_mimu")Reading layer `mmr_polbnda_adm2_250k_mimu' from data source
`C:\imranmi\imran's data sc\R-ex\R-Ex5\data\geospatial3' using driver `ESRI Shapefile'
Simple feature collection with 80 features and 7 fields
Geometry type: MULTIPOLYGON
Dimension: XY
Bounding box: xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS: WGS 84class(mmr_shp_mimu_2)[1] "sf" "data.frame"The Shape file for admin2 level map, is an SF object, with geometry type: Multipolygon.
st_geometry(mmr_shp_mimu_2)Geometry set for 80 features
Geometry type: MULTIPOLYGON
Dimension: XY
Bounding box: xmin: 92.1721 ymin: 9.696844 xmax: 101.17 ymax: 28.54554
Geodetic CRS: WGS 84
First 5 geometries:st_crs(mmr_shp_mimu_2)Coordinate Reference System:
User input: WGS 84
wkt:
GEOGCRS["WGS 84",
DATUM["World Geodetic System 1984",
ELLIPSOID["WGS 84",6378137,298.257223563,
LENGTHUNIT["metre",1]],
ID["EPSG",6326]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8901]],
CS[ellipsoidal,2],
AXIS["geodetic longitude",east,
ORDER[1],
ANGLEUNIT["degree",0.0174532925199433]],
AXIS["geodetic latitude",north,
ORDER[2],
ANGLEUNIT["degree",0.0174532925199433]]]Next, I will check the unique district names in this shape file (admin2)
unique_regions_mimu2 <- unique(mmr_shp_mimu_2$DT)
unique_regions_mimu2 [1] "Hinthada" "Labutta"
[3] "Maubin" "Myaungmya"
[5] "Pathein" "Pyapon"
[7] "Bago" "Taungoo"
[9] "Pyay" "Thayarwady"
[11] "Falam" "Hakha"
[13] "Matupi" "Mindat"
[15] "Bhamo" "Mohnyin"
[17] "Myitkyina" "Puta-O"
[19] "Bawlake" "Loikaw"
[21] "Hpa-An" "Hpapun"
[23] "Kawkareik" "Myawaddy"
[25] "Gangaw" "Magway"
[27] "Minbu" "Pakokku"
[29] "Thayet" "Kyaukse"
[31] "Maungdaw" "Mrauk-U"
[33] "Sittwe" "Thandwe"
[35] "Hkamti" "Kale"
[37] "Kanbalu" "Katha"
[39] "Kawlin" "Mawlaik"
[41] "Monywa" "Naga Self-Administered Zone"
[43] "Sagaing" "Shwebo"
[45] "Tamu" "Yinmarbin"
[47] "Kengtung" "Monghsat"
[49] "Tachileik" "Hopang"
[51] "Kokang Self-Administered Zone" "Kyaukme"
[53] "Lashio" "Matman"
[55] "Mongmit" "Muse"
[57] "Pa Laung Self-Administered Zone" "Danu Self-Administered Zone"
[59] "Langkho" "Loilen"
[61] "Pa-O Self-Administered Zone" "Taunggyi"
[63] "Dawei" "Kawthoung"
[65] "Mandalay" "Meiktila"
[67] "Myingyan" "Nyaung-U"
[69] "Pyinoolwin" "Yamethin"
[71] "Mawlamyine" "Thaton"
[73] "Det Khi Na" "Oke Ta Ra"
[75] "Kyaukpyu" "Myeik"
[77] "Yangon (East)" "Yangon (North)"
[79] "Yangon (South)" "Yangon (West)" There are 80 admin2 levels or districts in mmr_shp_mimu_2
Lets compare with our admin2 levels in our main dataset ACLED_MMR
unique_acled_regions2 <- unique(ACLED_MMR$admin2)
unique_acled_regions2 [1] "Maungdaw" "Bago"
[3] "Shwebo" "Kyaukme"
[5] "Pyinoolwin" "Muse"
[7] "Sittwe" "Yinmarbin"
[9] "Thaton" "Yangon-North"
[11] "Pa-O Self-Administered Zone" "Hpapun"
[13] "Kyaukpyu" "Yangon-West"
[15] "Mongmit" "Bhamo"
[17] "Mrauk-U" "Yangon-East"
[19] "Yangon-South" "Monywa"
[21] "Gangaw" "Pathein"
[23] "Katha" "Taungoo"
[25] "Kanbalu" "Lashio"
[27] "Mawlamyine" "Myitkyina"
[29] "Kawkareik" "Loilen"
[31] "Mandalay" "Kawlin"
[33] "Kyaukse" "Magway"
[35] "Meiktila" "Pakokku"
[37] "Taunggyi" "Tamu"
[39] "Nay Pyi Taw" "Mohnyin"
[41] "Kale" "Det Khi Na"
[43] "Myingyan" "Loikaw"
[45] "Matupi" "Pyay"
[47] "Sagaing" "Myeik"
[49] "Dawei" "Thayarwady"
[51] "Thandwe" "Mawlaik"
[53] "Bawlake" "Pyapon"
[55] "Hinthada" "Thayet"
[57] "Pa Laung Self-Administered Zone" "Mindat"
[59] "Hkamti" "Kokang Self-Administered Zone"
[61] "Hpa-An" "Danu Self-Administered Zone"
[63] "Myawaddy" "Maubin"
[65] "Hakha" "Falam"
[67] "Minbu" "Monghsat"
[69] "Puta-O" "Hopang"
[71] "Nyaung-U" "Kawthoung"
[73] "Yamethin" "Yangon"
[75] "Myaungmya" "Mong Pawk (Wa SAD)"
[77] "Oke Ta Ra" "Matman"
[79] "Kengtung" "Naga Self-Administered Zone"
[81] "Labutta" "Langkho"
[83] "Tachileik" I will write a simple function below to identify the discrepancies between the shape file and our state/district names in our main dataset.
# Find the unique region names that are in 'unique_acled_regions2' but not in 'unique_regions_mimu2'
mismatched_admin2 <- setdiff(unique_acled_regions2, unique_regions_mimu2)
if (length(mismatched_admin2) > 0) {
print("The following region names from 'acled_mmr' do not match any in 'mimu2':")
print(mismatched_admin2)
} else {
print("All unique region names in 'acled_mmr' match the unique region names in 'mimu2.'")
}[1] "The following region names from 'acled_mmr' do not match any in 'mimu2':"
[1] "Yangon-North" "Yangon-West" "Yangon-East"
[4] "Yangon-South" "Nay Pyi Taw" "Yangon"
[7] "Mong Pawk (Wa SAD)"Lets harmonize the names in both data files. I will re-save it to a new data set called ACLED_MMR_1
Fixing our admin 1 names.
ACLED_MMR_1 <- ACLED_MMR %>%
mutate(admin1 = case_when(
admin1 == "Bago-East" ~ "Bago (East)",
admin1 == "Bago-West" ~ "Bago (West)",
admin1 == "Shan-North" ~ "Shan (North)",
admin1 == "Shan-South" ~ "Shan (South)",
admin1 == "Shan-East" ~ "Shan (East)",
TRUE ~ as.character(admin1)
))Fixing our admin 2 names.
ACLED_MMR_1 <- ACLED_MMR_1 %>%
mutate(admin2 = case_when(
admin2 == "Yangon-East" ~ "Yangon (East)",
admin2 == "Yangon-West" ~ "Yangon (West)",
admin2 == "Yangon-North" ~ "Yangon (North)",
admin2 == "Yangon-South" ~ "Yangon (South)",
admin2 == "Mong Pawk (Wa SAD)" ~ "Tachileik",
admin2 == "Nay Pyi Taw" ~ "Det Khi Na",
admin2 == "Yangon" ~ "Yangon (West)",
TRUE ~ as.character(admin2)
))Checking if our changes are successful.
# Get unique admin 2 district names from 'ACLED_MMR_1'
unique_acled_regions2 <- unique(ACLED_MMR_1$admin2)
# Get unique district names from 'mmr_shp_mimu_2'
unique_map_regions_mimu2 <- unique(mmr_shp_mimu_2$DT)
# Find the unique district names that are in 'unique_acled_regions2' but not in 'unique_map_regions_mimu2'
mismatched_regions2 <- setdiff(unique_acled_regions2, unique_map_regions_mimu2)
if (length(mismatched_regions2) > 0) {
print("The following district names from 'acled_mmr_1' do not match any in 'mmr_shp_mimu_2':")
print(mismatched_regions2)
} else {
print("All unique district names in 'acled_mmr_1' match the unique district names in 'mmmr_shp_mimu_2.'")
}[1] "All unique district names in 'acled_mmr_1' match the unique district names in 'mmmr_shp_mimu_2.'"Lets do a sample plot to see how our country map looks like at the admin2 (districts) level.
plot(mmr_shp_mimu_2)
4.3.2 Data Wrangling
For the purposes of plotting choropleth maps, I will first create attributes subsets to summarise the number of incidents and fatalities, grouped by year, admin region, event type and sub event type.
Data2 <- ACLED_MMR_1 %>%
group_by(year, admin2, event_type, sub_event_type) %>%
summarise(Incidents = n(),
Fatalities = sum(fatalities, na.rm = TRUE)) %>%
ungroup()Checking the total sum of incidents and fatalities
total_incidents2 <- sum(Data2$Incidents)
total_fatalities2 <- sum(Data2$Fatalities)
total_incidents2[1] 57198total_fatalities2[1] 57593Next, I will do a spatial join between my shape files and attribute files
ACLED_MMR_admin2 <- left_join(mmr_shp_mimu_2, Data2,
by = c("DT" = "admin2"))Removing the variables I don’t require.
ACLED_MMR_admin2 <- ACLED_MMR_admin2 %>%
select(-OBJECTID, -ST, -ST_PCODE)class(ACLED_MMR_admin2)[1] "sf" "data.frame"Next, I will just double check that total sum of incidents and fatalities in our SF files are correct as per our original datasets.
total_incidents_check <- sum(ACLED_MMR_admin2$Incidents)
total_fatalities_check <- sum(ACLED_MMR_admin2$Fatalities)
total_incidents_check [1] 57198total_fatalities_check[1] 575934.4 Choropleth Maps
For the module on Cluster & Outlier Analysis (LISA), for the first tab, I will be plotting both Choropleth and Proportional Symbol Maps, as part of the initial Interactive Exploratory Analysis..
4.4.1 Choropleth map of Incidents & Fatalities by Admin2 level (by District)
The below codes will be used to create the choropleth maps.
tm_shape(mmr_shp_mimu_2) +
tm_borders() + # Draws borders for all regions
tm_fill(col = "white", alpha = 0.5, title = "Background") +
ACLED_MMR_admin2 %>%
filter(year == 2023, event_type == "Battles") %>%
tm_shape() +
tm_fill("Fatalities",
n = 5,
style = "quantile",
palette = "Reds") +
tm_borders(alpha = 0.5)
tm_shape(mmr_shp_mimu_2) +
tm_borders() + # Draws borders for all regions
tm_fill(col = "white", alpha = 0.5, title = "Background") +
ACLED_MMR_admin2 %>%
filter(year == 2021, event_type == "Violence against civilians") %>%
tm_shape() +
tm_fill("Incidents",
n = 5,
style = "equal",
palette = "Reds") +
tm_borders(alpha = 0.5)
Adding Interactivity by using tmap_leaftlet()
data_filtered <- ACLED_MMR_admin2 %>%
filter(year == 2022, event_type == "Battles")
tm_map <- tm_shape(mmr_shp_mimu_2) +
tm_borders() + # Draws borders for all regions
tm_fill(col = "white", alpha = 0.5, title = "Background") +
tm_shape(data_filtered) +
tm_fill(col = "Incidents", n = 5, style = "equal", palette = "Reds", title = "Incidents") +
tm_borders(alpha = 0.5)
tmap_leaflet(tm_map)